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Abstract 

The eccentric luminosity variation of quasars is still a mystery. Analytic results of this behavior 
ranged from multi-periodic behavior to a purely random process. Recently, we have used nonlinear 
time-series analysis to analyze the light curve of 3C 273 and found its eccentric behavior may 
be chaos [L. Liu, Chin. J. Astron. Astrophys. 6, 663 (2006)]. This result induces us to look 
for some nonlinear mechanism to explain the eccentric luminosity variation. In this paper, we 
propose a simple non-linear accretion disk model and find it shows a kind of chaotic behavior 
under some circumstances. Then we compute the outburst energy AF, defined as the difference of 
the maximum luminosity and the minimum luminosity, and the mean luminosity {F) . We find that 
AF ^ {F)°' in the chaotic domain, where a ~ 1. In this domain, we also find that (F) ~ M^'^, 
where M is the mass of central black hole. These results are confirmed by or compatible with some 
results from the observational data analysis [A. J. Pica and A. G. Smith, Astrophys. J. 272, 11 
(1983); M. Wold, M. S. Brotherton and Z. Shang, Mon. Not. R. Astron. Soc. 375, 989 (2007)]. 
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I. INTRODUCTION 

Since the discovery of quasar in 1963 the luminosity variation has played an important 
role in our understanding of its nature. Although it has been subjected to extensive analysis, 
there is no generally accepted method of extracting the information in the light curve. 
Results of analysis of the light curve ranged from multi-periodic behavior to 
a purely random process On the basis of these analysis, many theories which 

ranged from periodic mechanisms , llo| to superimposed random events such as the so- 
called Christmas tree model have been proposed. Then whatever does this seemingly 
random light curve tell us? Could this seeming randomness be some behavior other than 
multi-periodic or purely random? 

Also in 1963, Edward Lorenz published his monumental work entitled Deterministic Non- 
periodic Flow 12|. In this paper, he found a strange behavior which can appear in a de- 
terministic non-linear dissipative system, which seems random and unpredictable, and is 
called Chaos. Chaotic behavior is not multi-periodic because it has a continuous spectrum. 
Useful information can not always be extracted from the power spectrum of chaotic sig- 
nal. On the other hand chaotic behavior is not random either because it can appear in a 
completely deterministic system. The concept of attractor is often used when describing 
chaotic behaviors. As the dissipative system evolves in time, the trajectory in state space 

may head for some final region called attractor. The attractor may be an ordinary Euclidean 

( I 

object or a fractal [13] which has a non-integer dimension and often appears in the state 
space of a chaotic system. For many practical systems, we may not know in advance the 
required degrees of freedom and hence can not measure all the dynamic variables. How can 
we discern the nature of the attractor from the available experimental data? Packard et al. 
[l^ introduced a technique which can be used to reconstruct state-space attractor from the 
time series data of a single dynamical variable. Moreover, the correlation integral algorithm 

can be used to determine the 

dimension of the attractor embedded in the new state space. These techniques constitute a 
useful diagnostic method of chaos in practical systems. 

Therefore, the above-mentioned diagnostic methods have been used to analyze the light 
curve of 3C 273, and it has been found that the eccentric luminosity variation of 3C 273 
may be a kind of chaotic behavior [l^. This result tells us that non-linear may play a 
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very important role in the nature of quasar. Then whatever role does the non-linear play? 
Could the non-linear mechanisms help us understand the eccentric behavior in the luminosity 
variation of quasars? 

For the sake of keys, we in this paper propose a simple non-linear accretion disk (NAD) 



model, by borrowing basic ideas from the standard accretion disk (SAD) |17l. Il8[|. We then 
find that with some parameters the disk shows a kind of chaotic behavior. Although our 
model is based on the SAD, there are two main differences between them. First, SAD is 
a stable disk which can not be used to explain the eccentric luminosity variation. Second, 
SAD uses the complex differential equations of fluid dynamics to describe its evolution and 
NAD uses iterated equations. The latter simplifies the problem to a great degree. In the 
following sections, you will see very simple non-linear terms can produce extremely complex 
light curves. What's more, NAD is more natural than many other models on luminosity 
variation because our world is a non-linear world and the NAD aims to explore the non- 
linear effects on the quasar. Overall, our results provide a new view of the nature of the 
eccentric luminosity variation and may be helpful in the future study of quasars. 

II. DESCRIPTION OF THE CHAOTIC ACCRETION DISK MODEL 

There are always two kinds of factors which can be classified by their contrary contribu- 
tions to the behavior of the non-linear chaotic systems. For example, in the famous Logistic 
model 19|, one factor is related with the plenty food and/or space etc., and the other is 
related with the overpopulation and disease. The formal can make the population grow, and 
the latter can make the population decrease. The interaction of the two factors can lead 
to the chaotic variation of the population under some circumstances. The other example 



is the Lorenz model 121]. In this model, the two kinds of contrary factors are the heating 
of the earth, and the viscous stresses and the gravity. The heating of the earth can make 
the air rise. However the viscous stresses and the gravity always prevent the air rising. The 
interaction between them will lead to the chaotic behavior of the air. 

We note that the standard accretion disk (SAD) also has two kinds of contrary factors 
as in the above models. On the one hand, the viscous stresses will make the gas rotate 
slowly around the black hole. On the other hand, the gas in the disk will flow towards black 
hole, for its centrifugal force and gravity is not in equilibrium due to viscous stresses. The 
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gravitational potential energy will convert to the kinetic energy on the process of flowing 
towards the black hole which makes the gas rotate fast. Thus the viscous stresses and the 
gravity are the two contrary factors in the SAD. This is the start point of our non-linear 
accretion disk (NAD) model. By borrowing these basic ideas, we now construct the NAD 
model as follows. 

We consider an accretion disk which has only two layers of gas, as in Fig. [H It is supposed 
that (a) the thickness of each layer is very small that the tangential velocity in each layer 
is slightly different from its average tangential velocity. That is to say, we only need one 
velocity to describe the circular motion of one layer; (b) the height of disk H is very small 
that the motion is almost two-dimensional; (c) the average tangential velocity of outer layer 
is a constant and equals its average Kepler's velocity V . Due to the viscous stresses, the 
average tangential velocity of inner layer u does not equal its average Kepler's velocity f/, 
and is a time variable; (d) the inner layer interacts with outer layer by viscous stresses, but 
it does not interact with black hole; (e) during small interval Ar, the mass Am flows across 
any circle (all the centers of circles are on the center of black hole) is the same. If f/^ > m^, 
the centrifugal force is less than the gravity and the gas will flow towards black hole. Define 
Am > in such a case. If f/^ < m^, the gas will flow far from the black hole and Am < 0. 
When f/^ = M^, Am = 0. Thus, we simply suppose that Am oc U"^ — (f) there is a 
boundary layer between outer layer and inner layer. The thickness of this boundary layer 
h is very small comparing with the radius r of the inner layer's outer edge. When gas in 
the outer layer flows across the boundary layer, the average tangential velocity will change 
from V to U] (g) Any relativistic effect is not considered here. All above are the basic 
assumptions of our NAD model. By these assumptions, the equations of motion are derived 
as follows. 

At the n*^ moment, the viscous stresses acting on the inner layer is, 

/n = 27rrf//i — - — , 

where is the average tangential velocity of inner layer at the n*"^ moment and /x is the 
viscosity. Then we have, 

iTxrHjj, = m- 



h At 

where m is the mass of inner layer, Mn+i~ is the average tangential velocity of inner layer at 
(n + 1)*^ moment due to the action of viscous stresses, and Ar is the time interval between 



FIG. 1: The NAD model. The black dot at the centre of this figure represents a black hole and the 
shadow region is the boundary layer. Broad arrays represent the flow of mass towards black hole. 

n**^ and (n + 1)**^ moments. Thus it can be obtained that, 

27rrHfiAT 

Un+l =Un [Un-V). (1) 

When < f/^, the mass Am^ in inner layer will flow into black hole, which causes 
the same mass to flow across boundary layer into inner layer. When the gas flows across 
boundary layer, its velocity will change from V to U . According to the assumption (e), 



m V 

where C is a dimensionless parameter. At the (n + 1)*'^ moment, the average tangential 
velocity due to the flow of mass towards black hole is 

m m 

According to Eq. ([2]), we have 

When > f/^, there is only mass flowing out of inner layer, which will not change the 
average tangential velocity of inner layer. Thus, 

Mn+l"^ = Mn- (4) 

Finally, the average tangential velocity of inner layer at the (n + 1)**^ moment is 



According to Eqs. ([T]), ([3]) and (jlj), we have 



Un+l 



Mn + C(l - uJ/U^){U - M„) - {27irHfiAT/hm){u,, - V) when u^^ < 
Ujy — {27TrH fiAr / hm){un — V) when > U'^ 



Normahzing above equations, we have 

( u7, + 0(1- u;;^) {1- u;;) - A{u7, - B) when u;;^ <1 

Un+l = \ _ _ _ (5) 

I Mn — A{un — B) when > 1 

where = u^/U and dimensionless parameters 

We assume that the whole heat produced by viscous stresses was radiated out along 
the direction which is perpendicular to the accretion disk. Thus the luminosity at the n*^'^ 
moment is 

(u - V)'^ 

Fn = fn- {un -V) = 27rrg/i ^ " • 

h 

Normalizing above equation, we have 

F: = A{u;;-Bf, (7) 

where the dimensionless luminosity is 

mU"^ 

From Eq. ([6]), we can see that A is related to the viscous and geometry properties of the 
NAD. For avoiding the singularities in NAD model, the inequality 

dUn+l 



< 1 



must be satisfied when > 1. Then we can obtain that 

0<A<2. 

Under the assumption (a) and Eq. (Ill), we have 



^ _ /3r, I6GM 



c^r 
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where Tg is the Schwarzschild radius, c is the velocity of hght and M is the mass of the black 
hole. Clearly, B is related to the mass of the black hole and the geometry properties of the 
NAD and 

< < 1. 

The parameter C can also be related to some physical quantities. The distance that the 
mass Am at the outer edge of inner layer travels towards black hole during Ar is. 

As = Arl 

2r 



Then we can be obtained that 



A,n npHU^ Arl^^_^,^ 



m m 



Comparing Eqs. ([8]) with ([2]), we find that 



C 



We can see that C is only related to the geometry properties of the NAD. Here we assume 
that the mass of infiow or outfiow during At do not exceed the mass of the inner layer, that 
is, 

-1<^<1. (9) 

m 

For the inequality Amn/m < 1 is always satisfied when iT^ < 1, we choose 

< C < 1 

in the following discussion. From above discussions, we can see that only A is related to the 
viscous property of NAD, and only B is related to the mass of black hole. Thus by tuning 
corresponding parameter, we can easily know the effect of the variation of some physical 
property on the behavior of the NAD and its light curve. Some results are given in the 
following sections. 

III. RESULTS 

In this section, some properties of our model will be discussed by numerical computation. 
First, we will discuss the long-term chaotic behavior of NAD with different viscous properties. 
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FIG. 2: (a) and (b) are bifurcation diagrams of NAD model when B = 0.01 and C = 0.8. (a) is 
computed by picking a value of A and a starting point of = 1; iterating Eq. ([5]) 1000 times to 
allow the trajectory to approach the attr actor and plotting the next 5000 values of u. Then (b) is 
computed by Eq. (c) is the average tangential velocity nn of inner disk as a function of time 
and is computed by picking A = 1.95, B = 0.01 and C = 0.8, iterating Eq. ([5]) 1000 times from 
initial point uq = 1 and then plotting the next 100 values of u. For clearly showing its chaotic 
behavior, all the values of Un are minus 0.9 if they are greater than zero; conversely, they are plus 
0.9. (d) is the light curve computed by Eq. ([7]) also with A = 1.95, B = 0.01 and C = 0.8. 



As discussed in the above section, the viscous properties of NAD are just related to the 
parameters A. Thus in our following discussion B and C are fixed and A is an adjustable 
parameter. Fig. [2^ and Fig. [2)d are typical bifurcation diagrams when B = 0.01 and C = 0.8. 
They are very similar to the Logistic bifurcation diagram 20| . From these figures, we can see 
that the behavior of NAD will become more and more complex with the increase of viscosity. 
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FIG. 3: (e) and (f) are correlation integrals when A = 1.95, B = 0.01 and C = 0.8. They are 
computed by using 5000 data after 1000 iterations as in Fig. [2]^a)-(d), and the corresponding 
correlation dimension is given in each plot. 

Specially, there is a piece of vague area on the right side of each figure. The behavior of u in 
Fig. [2^ (or F in Fig. [Sjo) in this area will be chaotic. Fig. [2t is the evolution of the average 
tangential velocity u of the inner layer and Fig. 2d is the light curve generated by NAD 
model. Both of them are plotted when A = 1.95,5 = 0.01,C = 0.8, with which the system 
shows a kind of chaotic behavior. Here we must note that the choice of the parameters is 
not arbitrary. Any choice of the parameters should guarantee that Eq. (Q must be satisfied. 
By many numerical computations, we find that when A is greater than about 1.9 and C is 
greater than about 0.6, Eq. Q is satisfied and the luminosity is chaotic. 

The correlation integral [l5|] is also used to analyze our model. Let Xi,X2, . . . ,Xn be 
samples of a physical variable [u or F in our model) at the i*'^ moment, i = 1, 2, . . . , N. The 
correlation integral is defined as 

^ N N 

^°^(^) = NrN^^ ^ ^(R-|X,-X,|), 
^ ^ i=l j=lj^i 



where 9{x) is the Heaviside function, 

e{x) = 

In Ref. Il5|, they have studied many chaos models and found that 

Cor(R) oc R^, 



1 X > 
X < 0. 
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FIG. 4: The outburst energy AF as a function of the mean luminosity Fm- Here, we choose four 
groups of parameters. In each group, parameters A and C are fixed, whereas B varies in the interval 
(0,1), so long as the light curve is chaotic and Eq. ([9]) is satisfied. The samples of the light curve are 
produced as in Fig. [2] and the results are plot with lines and different symbols. The corresponding 
parameters and the slope of the line a are given as follows. Point, A = 1.95, C = 0.8, a ~ 0.93; cross, 
A = 1.97, C = 0.8, a ^ 0.96; plus, A = 1.98, C = 0.8, a ^ 0.97; star, A = 1.98, C = 0.6, a ^ 0.97. 

where D is called correlation dimension. Strictly D is not the dimension of the attractor, 
but is very close to it. According to chaos and non-linear dynamics, if an attractor for a 



dissipative system has a non-integer dimension, then the attractor is a chaotic attractor (20 1. 
Therefore, the correlation integral is a useful diagnostic tool of chaos. However, itself does 
not have any physical meaning. Fig. [3fe and Fig. [3]F are typical results of correlation integral 
when A = 1.95, B = 0.01, C = 0.8. From these results, we can see that attractors in the 
state space of u and F are indeed chaotic attractors. 

We then discuss the relationship of the outburst energy and the mean luminosity. The 
outburst energy AF is defined as the difference between maximum luminosity F(max) and 
minimum luminosity F(min), that is. 




AF = F(max) — F(min). 



This definition is the same as in Ref. l2l|. Here, we fix the parameters A and C and compute 
the luminosity with different B. For a guarantee of satisfying Eq. ([9]) and the chaotic 
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FIG. 5: The mean luminosity as a function of the parameter B. The parameters, the samples of 
the light curve and the corresponding symbols are the same as Fig. HI The line in the plot is 
(F) = A.3B + 2. 

behavior in the light curve, we choose A > 1.9 and C > 0.6, as was stated above. Then we 
compute the outburst energy and and the mean luminosity (F), and find that 



AF ~ {FY 



(10) 



where a ~ 1, as in Fig. |H This result is very similar with the observational facts found in 
the Fig. 7 of Ref. m 



where they found that for most sources it appears that the outburst 
energy scales with the mean luminosity. Additionally, we should note that B is related to 
the mass of the black hole. Thus this find suggests that the mass of black hole may not 
affect the nature of the properties of the luminosity variation. Quasars with big black holes 
will obey the same rules of the luminosity variation as the small ones. By fixing A and C, 
the relationship of mean luminosity and the mass of black hole is also discussed. We find 
that the mean luminosity scales with B, as in Fig. [51 With Eq. IQ, we can conclude that. 



(F) ~ M^/^. 

Combining Eq. ffTOl) with Eq. ffTTl) . we immediately obtain that. 



(11) 



(12) 
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where (3 ~ 1/2. Ref. [22 have reported that it is evident that the sources displaying largest 
variability amplitudes have, on average, higher black hole masses, although there is no linear 
relationship between them. Eq. (1121) is compatible with their results. 



IV. DISCUSSION AND CONCLUSION 

In this paper we propose a non-linear accretion disk model (NAD) which can be used to 



describe the chaotic behavior observed in the light curve of quasar 3C 273 16|. We note 
that the tangential velocity of accretion disk is influenced by two factors. One is the viscous 
stresses and the other is the flow of mass towards black hole. Under some circumstances, one 
of them may reduce the velocity and the other may increase it. Because of non-linear, the 
winner of gaining upper hand among them always vary with time. Then a kind of irregular 
oscillation of accretion disk would come out. It flnally leads to the chaotic luminosity 
variation. 

Although any complex multi-periodic mechanism and unnatural random event is not 
included, very simple non-linear terms in our model can also produce extremely complex 
behavior of the accretion disk. Meanwhile, for understanding the underlying non-linear na- 
ture of the chaotic behavior of the light curve, we ignore many physical details, including 
some may be very irnportant to the behavior of the quasars. For example, the effects of 



the thermodynamics 23l |. the interaction between the inner edge of the inner layer and the 



black hole [2j], and the details of the dynamics happening in the boundary layer are not 
considered here. Additionally, we use the iterated equations instead of the complex equa- 
tions of fluid dynamics. All above-mentioned facts simplify the problem greatly. However, 
the simpliflcation does not weaken the validity of our model on explaining the eccentric 
luminosity variation. The model can shows two rules in the chaotic domain. One is the 
outburst energy AF ~ {F)"", where {F) is the mean luminosity and a 1; the other is 



{F) ~ M°'^, where M is the mass of central black 



ho 



compatible with the observational data analysis 



21 



22|. 



e. There rules are conflrmed by or 



The standard accretion disk model (SAD) is always be recognized as a standard picture 



of the AGNs, for it can produce extreme energy that be observed in the AGNs 25|. Our 
NAD model just borrows the basic idea of the SAD. Thus, we conclude that the chaotic 
behavior presenting in NAD model can be seen in many eccentric light curves of the AGNs. 
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Some authors have reported that other AGNs may have chaotic behavior in their eccentric 
hght curves 26^. However, we hope that there will be more evidences to support or oppose 
the assertions of our model, especially the quantitative rules predicted in this paper. 

Overall, our results reveal that the non-linear plays a key role in the cause of the eccen- 
tric luminosity variation of quasars. The omnipresent non-linear is very important in our 
understanding of this complex world, as chaos and non-linear dynamics tells us. Thus it is 
nature to study the luminosity variation of quasars from the view of non-linear. And we 
hope chaos and non-linear dynamics may be helpful in the further study of quasars. 
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